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Abstract It is well known that mode coupling theory (MCT) leads to a two-step power-law time decay in dense 
simple fluids. We show that much of the mathematical machinery used in the MCT analysis can be taken over to 
the analysis of the systematic theory developed in the Fundamental Theory of Statistical Particle Dynamics [1]. We 
show how the power-law exponents can be computed in the second-order approximation where we treat hard-sphere 
fluids with statics described by the Percus-Yevick solution. 



1 Introduction 

The long time dynamics of dense classical fluids is complicated. There is experimental evidence of power-law 
decays as a function of time. One approach used to describe this evolution is Mode-Coupling Theory (MCT). MCT 
has been an extremely controversial theory. The problem is that there is no systematic derivation. MCT has been 
cobbled together piece by piece leading to its present form. There has been no systematic way of deriving, correcting, 
or extending the basic description. 

Despite these criticisms MCT has led to a long time scenario which seems to match a lot of the physical opera- 
tions in a dense fluid. It has been argued by |2] that MCT is the "mean field" theory goveming glassy behavior and 
the basic results are more general than found in MCT. 

The most contentious issue in judging MCT is the existence of an ergodic-nonergodic (ENE) transition. It is 
commonly believed that there are no ENE transitions in the classes of models studied. However there is almost a 
transition. In the ideal MCT scenario there is an ENE transition at some critical packing fraction T]*. MCT implies 
a rather elaborate slow dynamics on either side of the critical density. There is, for 77 < T]*, a three step process |i3J. 
The equilibrium density autocorrelation function first exhibits power-law decay 

Gppiq,t) = (p_g(o)pg(f))c'-/g+ v"^ 

followed by the von Schweidler decay 

(p_,(0)p,(f))c-^/,-V> 
and finally exponential (possibly stretched) decay 

{p.,{0)p,{t)),^e-('M'. 
A full analysis of the theory gives expressions for the exponents 

r{i-af ^ ni+b)^ 

r(l-2a) r{l+2b)' 

where A is a material dependent quantity. 
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Fig. 1 Generic one loop contribution to Ssb (1,2). 

In FTSPOfT] and SDENE|4] one of us introduced a new theory for the dynamics of collections of classical 
particles. It is shown in |4] that at second order in perturbation theory one is led to a description which has similarities 
with MCT. 

In this paper we take the analysis of the kinetics beyond the analysis of the location of the transition. We explore 
the slow kinetics for the regimes bellow and above the transition. It was shown in |3, 0] that the density-density 
cumulant Gpp {q,t) satisfies the kinetic equation 



dG 



PP 



dt 



{q,t) = -Dq^pS{q)-^Gppiq,t) + Dq^ ds p^pZBBiq, 



t-s)Gpp{q,s 



(1.1) 



where j3 = 1 /hgT is the inverse temperature, p is the average density, D is the diffusion coefficient, S{q) is the static 
structure factor, and Ebb is the B-B component of the self-energy matrix for the theory of the two fields p and B 
where S is a response field. The important point is that Lbb can be determined systematically in perturbation theory 
and, at second order in the potential, is given by 



P^p'LBBiqut) = y dt' d{qi -k2- h)d(t - t' 



d_d_ 
dtd? 



Gpp {k2,t)Gpp (^3 ,t')+Gpp {k2,t)Gpp {ki,t' 



where K12 = {qi ■ ^2) /q\k\, = {q\ ■ k^)/q\k^, Gpp is the non-interacting density autocorrelation function and Gpp 
and Gpp are dressed propagators Recently it has been suggested |6] that it is important to include information in 
addition to the static structure factor to make accurate quantitative predictions about glassy dynamics. In principle, 
a systematic extension of FTSPD to higher orders in the effective potential can take such data into account. 

There are a number of differences between the MCT kernel and Lbb although the overall structure is similar. 
In the FTSPD the self-energy has the basic structure represented by the one-loop diagrams shown in Figure 1 . The 
triangles at the ends of the graph are the lowest order vertices. We study here the cases where the vertex functions 
are given by the non-interacting values. These vertices have a frequency dependence. Between the vertices we have 

several sets of products of matrix propogators. An interesting feature is that all of these propagators g|°', G,,, Gij, 
Gij satisfy similar fluctuation-dissipation theorems. For a full discussion see |0]- 

It has been understood for some time that there are corrections to MCT that consist of higher numbers of loops. 
The two loop correction is of the general form 




where now one has 4-point vertices and a product of these dressed matrix propagators. 

Our work here is to explore the long time solutions of the kinetic equations. Near an ENE transition we can 
assume certain analytic forms for the time dependence of Gpp{q,t). For MCT a rather elaborate mathematical 
machinery has been developed for performing this asymptotic analysis |8]. While none of this is new for MCT, we 
find it useful to present these procedures in some detail. The results of an analysis for MCT are well known; for the 
case of the FTSPD the results are all new and comparable to the MCT. We begin with some background information 
on FTSPD and MCT. We follow this by a review of the asymptotic MCT analysis. Finally, we apply the analyis to 
FTSPD and numerically determine the power-law and von Schweidler exponents for hard spheres using the Percus- 
Yevick approximation for the static structure factor. We compare these results to those of the corresponding mode 
coupling theory. 
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2 Reviews of FTSPD and MCT 

In this section we shall briefly review the origins of FTSPD and MCT and how one can use them to construct kinetic 
equations like dl.lb for correlation functions of interest. 

FTSPD 

The Fundamental Theory of Statistical Particle Dynamics is a framework from which mode-coupling-like kinetic 
equations can be derived, including systematic corrections to them. In Smoluchowski dynamics, FTSPD is devel- 
oped as a field theory whose fundamental field is the position variable, Starting with the equations of motion, 

Riit)=DFi{R) + 7]iit) 

together with the variance of the Gaussian noise term 

{r]i{t)r]j{t'))=2kBTDdijd{t~t'), 

we obtain a generating functional for a path 

rR{t)-H{t)c 

•Itn 



where R H is shorthand for 

poo 

and the integrated-over initial conditions are 



Ri{to)=RT 

with probability distribution given by Po(^,-''')- We can obtain the Martin-Siggia-Rose action for this system by 
exploiting the functional delta-function identity 

j [d7]i] 5 {Ri - DFi -rii)=J [dRi] detM; 5 (/?,■ - DFi - t],) = 1 , 

where detM, is the Jacobian 

detM, =det4^. 

oRi 

Inserting this identity into the generating functional and using the functional Fourier representation of the delta- 
function we obtain 



Z{H) = I [dm]d^Rf^P^{Rf^)^^^^-^^ j\fdt +RrHi 

= j [dr]i][dRi]d'Rf^Po{Rf^) detM,-5(/?,- -Z)F, - T],)exp Tif^/f +/?, .//, | 

= j[dni] [dR,] [dR^d'RfPoiRf) detM, exp ^,(/?, - DFi - 77,)^/^ | exp | -^J^ 7]fdt +Ri ■ 
We now combine the arguments of the two exponentials into a single integral: 

Z{H) = J [drii] [dRi] [dRt]d^Rf^PoiRf^) detM,- exp | jj' (^/^,(/?,- - DFi - r?,) - 5^^) dt + Ri ■ //, | 

= J [dr]i] [dRi] [dRi]d'Rf^PoiRf'^) detM,- exp | jf" (/^,-(/?,- - DFi) - '^^-^J£^§^) ' 
= j [drii][dRi][dRi]d^Rf^Po(R\°'^)dctMi 

^c^^^^j^^l^R,[R,^DFi)-^^^i±^^^-UBTDR'^ dt+Ri-H^ . 
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Now all other variables are independent of 77,, so we can perform the Gaussian integral and absorb the constant into 
the measure. We obtain 

Z{H) = j [dRi][dRi\d^Rf''Po{Rf^)Ae\.MiCj.x>^- j RikBTDRi - iRi{Ri ~ DFi)dt + Ri ■ Hi 

We use this to build the generating functional for all the paths: 

Zn (H) = J Y[[dRi] [dRi]d^RrPo {Rr )expl-J dt £ [RiksTDRi - iRi{Ri - DFi) -A)] + ^-H\, 

where we have rewritten the Jacobians as A'j = IndetM, and introduced the symbol (/) to represent the fields (con- 
structed from R, R) that we wish to couple to H. A more complete derivation with a treatment of the Jacobian is 
available in lU]. From the generating functional of the paths we obtain the grand canonical partition function 

N=0 ■ 

from which we construct the generator of cumulants 

W[H]=logZT[H]. 
The source fields Hp and Hg are coupled to the density 

p(l)=£5(xi-/?,(fi)), 

and response field 

B{1) = Z)£ {^,Vi + 0(O)V2} 5ixi -Riiti)) , 

i 

respectively. The value of 6 (0) is chosen based upon the discretization scheme we adopt for the path integral. The 
standard Stratonovich discretizatiorQ results in 0(0) = 1/2. We generate cumulants by functional differentiation; 



G,.,.....(1,...,£,£+1,...,„) = ^...^^...^H^[//] 
From the Legendre transform of W [H] one generates the irreducible vertices 



Hp=HB=0 



5 5 5 

The two-point vertices satisfy the Dyson equation 

k 

Expanding the vertex to second order in the potential we have i7i. = y^j^"^ + — Lij^. Plugging into the B-p compo- 
nent of Dyson's equation we construct the kinetic equation 

r)G ft 



_ ft 

^{q,t) = -Dq^pSiq)-'Gpp{q,t)+Dq^ ds pLBB{q,t ~ s)Gpp{q,s) , 



dt 

where D = ksTD and Gpp {q,0) = S{q). It is shown in [4] that approximations for Ebb can be generated as a power 
series in an effective interaction potential. 

' Using the Ito prescription is equivalent to making 0(0) = 0, resulting in a trivial Jacobian. The Stratonovich prescription is usually 
preferable since the ordinary chain rule of calculus remains unmodified when dealing with multiplicative noise. Either way, in FTSPD 
the Fokker-Planck equation and all coirelation functions resulting from the path integral are independent of the choice of discretization 
since the noise is not multiplicative. 
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MCT 

The most common path to Mode Coupling Theory is through the projection operator formalism of Zwanzig and 
Mori. We briefly summarize the development given in [9]. We start from the fundamental relationship between a 
function of phase-space variables A (r,/)) and the Hamiltonian M' = L//'^/2w+ jUi^jVifij)- 

dA 

— ^{A,J^} = i^A{t), (2.1) 
at 

where we have used the Poisson bracket 

' dA dB dA dB 



{A,B} = y. ^ ^ ^ ^ 

This equation defines the Liouville operator 

i m dri f^. dri dpi 
It follows that we may formally express the time evolution of A as 

A{t) = e'^'A{0) . 

The next step in the analysis is to identify A with a set of functions that we are interested in tracking. These are 
referred to as "slow variables." For the derivation of the basic MCT equation these are the density fluctuations and 
longitudinal current: A{t) = {5pk{t),j]^{t)). We also define an inner product 

(A,S) = (A^'B), 

which allows us to define the projection operator0 

3^=A{A,...){A,A)-K 

Note that A = A(0) above. The projections perpendicular to the slow variables are referred to as "fast variables." 
Separating the left-hand side of (12.1b into components parallel and perpendicular to A(0) allows us to express the 
time evolution as 

^ = iQA{t)+ f e'^'^'-'H^^f{T)dT + f{t), 
at Jo 

where 

in = (A,;^A)(A,A)-i 

and 

f{t) = e'^'-^^^'iil-I^)^A. 
Exploiting the orthogonality between /(/) and A(0) allows us to write 



dA 
dt 



= iQA - [' K{t-T)A{T)dT + f{t), 
Jo 

where the kernel is given by 

^(0 = (/(0),/(f))(A(0),A(0))-^ 
Introducing the correlation matrix C(f) = (A"'"(0)A(f)) = (A(0),A(f)) we have 



dC r' 

^=iQC{t)- / dTK{t-T)C(T). 
dt Jo 



^ The projection operator retains the same form when A is generalized to a set of functions. In this case, {A, A) refers to the matrix 
Aij = iAi,Aj) and = (A,A)-',/Aj,B)A,. 
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Specializing to our choice of slow variables, A = {5pk{t),j\{t)), we obtain 
from the p-j component of our correlation matrix, where 

Gpp{q,t) = {5p-q{0)5pq{t)) 

is the density-density cumulantEl and 



at mS{q) 

So far all of these manipulations have been exact. However, computing with ( 12.2b is not practical without simplifying 
the kernel K{t) = {R-q{0)Rq{t)). This is the point where the uncontrolled approximations of MCT are introduced. 
They consist of: 

- Replacing the time evolution operator ^je'"^' where ^2 projects onto the dominant "slow" 
mode of Rq. 

- Gaussian and convolution approximations for the correlation functions which allow them to be factorized into 
simpler correlation functions proportional to the static structure factor. 

In contrast to the FTSPD, it is not clear how one could go about systematically improving the MCT approximations 
for the memory function kemel. 



3 Review of MCT Asymptotic Analysis 

These approximations define mode coupling theory and lead to the following equatiorQ for the density autocorrela- 
tion function, Gpp{q,t): 



^Gpp(q,t)+Q\q)Gpp{q,t)+Q\q)S{q)l^ dt'M{q,t~t') — Gpp(q/}=0, 



where 



q^keT 



mS{q) 

M{q.t) = y ■0^V^{q,k)Gpp{k,t)Gpp{\q~klt) 

y2(q,k) = 4 • ^c{k)+q- (q - k)c(|q- k|)}2 
1 

pc{k) = 1 - ^ 



S{k) 

The important point is that M is a non-linear function of the density autocorrelation function. V{q, k) is the vertex in 
the theory. We may rewrite the equation in terms of the normalized density autocorrelation function ^{k^t) which is 
related to the autocorrelation function by Gpp {q,t) = S{q)^{q,t). The kinetic equation can then be written as 

^(^{q,t)+n\q)Hq,t)+n^{q) dt'M'{q,t-t')-^Hq,t')=0 . (3.1) 

^ It should be noted that the convention used in this section (and most of the MCT literature) for the cumulants is different 
from the convention we have employed in the discussion of FTSPD: G^J(g,?) = lim V"' (5p_,(0)5p,(f)) while G'^^'^{q,t) = 
lim A'^' (5p_f^(0)5p^(«)), where V and A' are the volume and number of particles, respectively. Accordingly the static structure factors 
are related in the same way: S'^'^iq) = pS'^^^iq). We continue to use the MCT convention when we discuss the MCT kinetic equation 
below, and switch back to the original convention when we discuss the FTSPD kinetic equation. 

We are presenting the Newtonian MCT equation. For Smoluchowski dynamics only one time derivative appears. 
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The corresponding kernel M'{q,t) is given by 

1 r d\ 



M'{q,t) = -j y2(q,k)(/)(^,f)0(|q-k|,O 



{2n) 

y2(q,k) = ^S{q)S{k)S{\vi ~V\){q- ^c{k)+q- (q - k)c(|q - . 

q 

The only relevant parameter for hard spheres is the number density, p. We are interested in extracting the long-time 
slow dynamics for this model. In the following sections all numerical calculations and plots are found using the 
Percus-Yevick approximation for the hard sphere sturcture factor. It is useful to work using the Laplace transform^ 
which diagonalizes ( 13. lb This leads to 

Z^^{q,z)-Z + ^^^{q,z)+^^^[M{q,t)]{z^{q,^) - 1) = , 
where we have applied the initial conditions ^{q,0) = 1 and (j){q,0) = 0. Rearranging terms we have 

^^^^^ = Q^^[Miq,t)]+z. (3.2) 
l-z(l){q,z) 



3.1 Heuristic Analysis of MCT: Time Dependence at Transition 

In the long-time regime we assume ^{q,z) — > oo as |z| — > 0. In this regime (|z| ^ 1) we may ignore the z on the right 
hand side of ( 13.2b and thereby cancel the factors of ^2. This gives u^ 

•^(^'^^ =^[M{q,t)]. (3.3) 



1 -z(/)(^,z) 



Next we apply the ansatz <^{q,t) = f{q) + (1 — f{q))^Wq{t) where \lfq{t) is assumed small and f{q) is a non- 
ergodicity parameter. The Laplace transform of the ansatz is given by 

</>(^,z) = — + (!-/(?)) V.W- 

z 

We now plug this into ( 13.3b and expand the left hand side to second order in ^fq{z): 

(7 + (1 -/.) V^.W - Ufwliz)) = ^[M{q,t)] . 
To apply the ansatz to the right hand side note that 

(/.(fc,0</>(|q-k|,0 = {h + {\-fkfWk{t)}{fq-k + {\-fq-kf%-k{t)} 

= fkfq-k+fkil - fq-k?Wq{t) + fq-kil ' fkfWkit) + (1 - ( 1 (f) V/,-,(f) . 

The two middle terms above may be combined into a single term in the integral. We now employ the notation of 
to make the following more compact. Let 

pS{q) f (Pk 



5p 



Cqk = ^^{l-fk)' (3.4) 

Ofk 

2 dfkdfp 



^ Our conventions are defined in Appendix B. 

^ Note that in MCT both Newtonian dynamics and Smoluchowski dynamics satisfy the same asymptotic equation. 
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We obtain 



^ ^ —-F,[fk] ] + {w,{z)-C,kWk{z)) +z{\~f,)Wq{z)~CqkpJ^[xifk{t)Wp{t)] = 



z v-u 



(3.5) 



where it is understood that we are integrating over repeated indices: Cq]^/]^ = J t/^kC^j^/t/(27r)^, etc. This expansion 
yields three conditions for the ideal metastable state. The first one, given by 



= F,[fk] 



determines fq. For hard spheres, this equation can be solved by iteration. A plot of fq using the Percus-Yevick 
approximation for the hard-sphere structure factor is available in Figure |2] The second and third conditions can 




ka 



Fig. 2 This plot shows ./j,^*"^, the non-ergodicity parameter, at the critical density tJmct = 0.515. 



be found by considering a further expansion of iffqiz) = v4^^ + v4^^' '^^^^ v4^Vv4^^ ^ as z —> 0. We plug the 
expansion into our equation and match orders. To the order of ^fq^\ we obtain 



v^'\z) = Cq,xif;^"{z), 



(3.6) 



where Cqk is defined by (13.4b . This equation determines the spatial dependence of i//^'' (z). We can solve this equation 
for y/q^^ (z) by finding an eigenvector of the operator 



Cqk — 



^y2(q,k)/,_,(l-/,)2 



with unit eigenvalue. This eigenvector is found at the critical density 111211 . Below the critical density (13.6b cannot 
support a nontrivial solution. Above the critical density the eigenvalue moduli are all less than one. Denote the right 
eigenvector associated with the unit eigenvalue by ej^, and the corresponding left eigenvector by e^. Furthermore, 
impose the conditions Y^k^k^k = 1 and Y.k^k^k^k{i — fk) = With this notation we have 



V4 =Aeq(j)v{z) , 
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Fig. 3 This plot shows e^'~^ , the right eigenvector of C^*"^, at the critical density t^mct = 0.515. The eigenvector is normalized 
according to the convention Y,k^k^k = 1- 

with A a thus far undetermined constant and (j)v{z) a function of z only. Plots of e^. and for the hard sphere system 
are available in Figures |3] and |4] respectively. Next, we match the second order terms in (13. 5t . We obtain 

(V-c,,)</>f = -(1 -/,Mv/^'']2+Q,,^[va1'V^'^] ■ 

Now multiply both sides by the left eigenvector, iq. We obtain 

= -e,(l -/,)z[v/^')]2 + e^Q,,if[v/l'Vp''] 

where we have used the condition Y^k^k^k^ki^ ^ fk) = ^ and A = eqCqkpCkSp in the third line. This gives us an 
equation for the z-dependent part of v'^^ : 

z(/»2(z)=A^[</>v'(0]- (3.7) 
Usin^ ^v{t) = t~", which implies ^v{z) = we obtain a quantitative expression for the power law index 

r^{\-a) ^^ r{\-2a) 

This implies 

r^{\-a) _ ^ _ , 

r{\-2a) -^-'cfqkpekep, 

where X has to be determined numerically. This analysis only applies near the transition density. The assumption of 
y/- = y/C) -|- y/(2) and matching the orders might fail in other circumstances. A question remaining in this analysis is 
the nature of the separation of ^fq{z) into first and second order terms. 




' The amplitude of the power law is absorbed in the parameter A. 
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Fig. 4 This plot shows e^^"^, the left eigenvector of C^^^, at the critical density tJmct = 0.515. The eigenvector is normalized according 
to the convention Y^k^k^k = 1- 



3.2 Ergodic-Nonergodic Transition 

Having demonstrated the approach corresponding to slow power-law decay, we can push harder for a fuller descrip- 
tion at the transition as well as above and below. We take first the case of working at the transition. To see if the 
mode-coupling approach supports an ergodic-nonergodic transition, we write (j>{q,t) = fq + gq{t), with gq{t) as 
f — > oo. This implies (l>{q,z) = fq/z + gq{z), withg^(z) less singular than 1 /z. Let us make this substitution in ( 13.2b : 



n\fq/z + gq{z)) 



l-z(fq/z+gq{z)) 

Ignoring terms less singular than 1 /z in the long time limit (z — > 0) we obtain 



1-/. 



Fqifk). (3.8) 



Therefore, we find that there is a non-ergodic state when the above equation supports a nontrivial solution for 
Solving this equation by iteration reveals that there is indeed an ergodic/non-ergodic transition for hard spheres at 
7] w 0.52, using the Percus-Yevick approximation for the static structure factor. 

Next let us determine how the non-ergodicity parameter changes when we perturb the density: p — > p + 5p. 
Write fq = fq +{i - fqYgq, where gq is "small", and substitute into the left hand side of ( 13.81) : 



U _ f,+{^-f,?gq _ fq 



The right hand side of ( 13.8b expands to 

op Ofk 2 dfkOfp 

where 5fk = (1 — fkYgk and Fq is evaluated at Using notation introduced in Section lTTI we obtain 

pp+5p ^fP+5p )^pP^Cqe + Cq,gu + Cqkpgkgp + 0{g\gp,p'). 



'"(/r'") = FP + ^5p + ^dfk +i4J^8f,df„ + o[g\gpy) , 
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Combining both sides we obtain an equation for g^: 

^9 + (1 - f^)g\ = Cqe + Cqkgk + Cqkpgkgp + 0{g\gp,p'^) . 
Rewriting in a more suggestive form, 

{^k-Cqk)gk = -{^-fq)g]+Cqe + Cqkpgkgp + 0{g\gP,P^) (3.9) 

At the transition density, p = Pc, Cqk has unit eigenvalue. Let us expand gi^ = +^[^' ^nd assume that is 

the unit eigenvalue of Cqk- In that case, we can write = gei^, with g a thus far undetermined small parameter. 
Plugging this form into (|3.9l l allows us to relate g and the separation parameter e = (p — pc)/pc- Expanding to 
second order in g, we obtain 

{8qk Cq,)gf^ = -{^-fP){g^q^Y+Cqe+Cqkpg^S^P + 0{g\gP,p') . 

(2) 

Multiplying both sides by the left eigenvector of Cqk to eliminate g^ ' yields 

= -g^eg{l-fP)e^^ + egCqS + g^e^Cqkpekep + 0{g\gp,p^) . 
Recalling the convention for the eigenvectors, Cg(l — fq )eq = 1, and the definition of /I = eqCqkpSkep, we find 

2 _ ^tjCqS 

^ 1-A ■ 

We also know that Cq is positive, since to leading order in the separation parameter Fq^^^ = (1 + £)Fq . Therefore 
our assumptions are consistent near the transition density as long as A < 1 and the small parameter in scales as 
\/e. Putting all of this together, the non-ergodicity parameter near the transition density can be written as 

U=fq+g{'^-Ufeq- 

A plot demonstrating the proper scaling behavior is shown in Figure|5] A similar analysis implies that for the largest 
eigenvalue Eq of Cqk, assuming the eigenvalue is one at the critical density, Eq = \ — 2g{\ — A), agreeing with the 
statements in lfl3ll . A plot of the maximum eigenvalue of Cqk as a function of the packing fraction is shown in Figure 
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3.3 Time regimes 

In this section we will analyze the asymptotic long-time behavior of (j){q,t). First we need to find out what happens 
at the critical density so we can match the other solutions to it in the limit e — 5- 0. The analysis of (j>{q,t) at the critical 
density is essentially the same as in Section [ITl but we include it here for completeness. In the non-ergodic (glass) 
case we have only one long time to consider. In the ergodic (liquid) case there are two time regimes: the power law 
decay as in the solid case followed by the von Schweidler decay. To analyze these regimes we rescale the time t by 
introducing a large constant (Oc such that 1 = (Oct. In this rescaling the power law decay corresponds to T ^ 1 and 
the von Schweidler decay to t » 1. 

3.3.1 Critical Density 

Here we are interested in the asymptotic behavior of ^{q,t) at large times. First we consider what happens at the 
critical density. As in Section lTTI we apply the ansatz (^{q,t) = f{q) + (1 — f{'i)YWq{^)- We obtain 

^ -^-?[/^]) + {^q{z)-Cqk^fk{z)) +z{\-fq)^fl{z)-CqkpJ^Wk{tWp{t)]=0, (3.10) 

where the first term (proportional to 1 /z) vanishes by definition. Again, considering a further expansion of ^iq{z) = 

y/q^^ + yi'^^ with yilp /yi'i^'^ — ;> as z — ;> 0, we plug into our equation and match orders. To the order of v4'^' 
obtain 

wl^\z)=Cqkw]^\z). 
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n 

Fig. 5 This plot shows f^^{ka =1.0), the non-ergodicity parameter at wavenumber ka = 1.0, as a function of the packing fraction. 
The dashed line shows a visual fit of A-^e + B to the data. This confirms our calculation of g in Section [l!2l 



Therefore, we may write Xffq (z) = Aeq^viz)- Plugging this back into ( 13.101) and multiplying by the left eigenvector 

)r (/)v(z): 



(2) 

of Cqic to eliminate Xffq , we obtain an equation for (/)v (z) 



3.3.2 Glass Side 

Next let us consider what happens when we are slightly away from the critical density on the glass side. Here we use 
the ansatz ^{q,z) = fq/z + (1 — fq)^gWq{z), with ^fq{t) — > as f — > oo. {g was computed in Section [l!2l ) Starting 
from (13.2b and expanding to second order in g. 



\ -^^[A] -Q^j +g{xifc,{z)-Cq,Xifk{z) - dCqkWkiz)) +g' {z{\~fq)wl{z)-CqkpJ^[Xifk{t)Xif,{t)])=0. 

According to our assumptions the first two terms in the 1 /z factor cancel. We are left with 

-\ +8ijc,{z)-CqU\lfk{z)-5Cqk\ifk{z))+g' {z{\-fq)\ifl{z)'Cqk„J^[Xifkmp{t)])=0. 

Let us re-write 5p /p as g^B. (B can be determined from the expression for g^ in Section 4. Note that e = 5p/p. 
Also note that dCqk = gDqk since 5Cqk scales as \f£.) We obtain 

g{xifq{z)-Cq,Xif,{z)) +g' (^z{\-fq)^fl{z)~Cq,pJ^[^fkmp{t)]-Dqm{z) - ^) = . 

Using the same method as before, we split i//^ = vl'' +8Wk''' ■ We find 

=AekWv{z)- 

Using this result and eliminating y/'^) by multiplying by the left eigenvector we obtain an equation for ^fv{z)'- 

zwl{z)-:k^Wl{t)]-D'^fv{.z) - - = , (3.11) 

z 
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1.2 




0.504 0.506 0.508 0.51 0.512 0.514 0.516 0.518 0.52 

Packing Fraction 

Fig. 6 This plot shows the maximum eigenvalue of C^^^ as a function of density. The lines are drawn at 30, 60, 120, 240, and 270 
iterations. We see the lines converging as the number of iterations increases. 



which is slightly different from (13.7b . Now change frequency variables to ^ =z/(Oc, with e — > 0, (Oc/s — > z — > 
and 1^ = z/(Oc constant. This implies .if [/(?)] (z) — )■ 0)~^ ^[f{T)]{Q, where T = (Oct. Rewriting according to this 
change of variables: 

We see that we obtain the same equation as before, but with z — > i^. Now let us make the ansatz V'v(i') = 1 + T"". 
This implies Yv{Q = l/^ +r{l - a)/!^^'". Substituting this in dTni l. 

,f r{i-af , , r(i-a) , 1 ^ ^ f r(i-2ay , , r(l-a) , l\ ^,f r{l-a) , 1^ B' 
In the limit ^ > 1, T < 1, we find that 

r(]_a)2 
r(i-2fl) ' 

provided < a < 1 /2. In the limit ^ 1 the solution goes to its final asymptotic value. Plugging in above allows 

us to obtain an equation for A in ^fj.^\t) = Aei^(j)v{t), removing some ambiguity from our solution. Evaluating the 
limit as e — )■ and matching the solution with what we found at the critical density implies that 

limg{l-f,fe,{l + T~") = {\-f,fAe,t-\ 

Since g e'/^, we must have 

gi/2^-a ^ constant x t~" . 

From the definition t = (Oct, 

£^/^(0c" = constant. 

From this we find that (Oc scales as e^/^". 
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3.3.3 Liquid Side 

Finally we consider what happens when we are slightly below the critical density (liquid side). Here we need to 
assume that a slight change from the critical density results in a long plateau region which is close to Following 
the same procedure as before we obtain a slightly different equation: 

^yUQ ^^[vlmQ +D'wv{Q + 1 = 0. 

Note the change in sign of the final two terms. This is due to the separation parameter being negative on the liquid 
side of the transition, implying g is negative there as well. However, the analysis for ^ 1 is basically the same as 
on the glass side, and we obtain 

r(]_-«)2 

r(l-2a) ■ 

This solution is valid as long as fo ^ f <C co~^, where to is a timescale governing the microscopic dynamics and (Oc 
was found in the previous section to be Qq x e'/^" with Qq constant. 

For 1^ <C 1 we consider the ansatz Yv{'^) = 1 — f*, which implies (//y(i^) = l/i^— F(l Plugging this 

into our equation, 

rf, . r(i+j.) r{i+bf \ ( r{\+b) r{\ + 2b) \ f r(i+fe) \ b' 

Taking the limit, and considering only the most singular terms, we obtain 

r(i+^ 
r(i+2^) • 

provided ^ > 0. The validity of our solution is based on the assumption that \z'^{q,z) \ ^ 1. We used this assumption 
when we expanded the fraction in (13.2b as a geometric series. Since our solution is singular as z — S> 0, this inequality 
will eventually be violated rendering our solution invalid. This will happen when 

where y = 1 /2a + 1 /2b. Therefore this von Schweidler part of the solution is valid when 

Thus we have analytically found the famous two step process of a power law to a plateau followed by the von 
Schweidler decay. 



4 Application to FTSPD 

It tums out that the asymptotic analysis of the FTSPD is remarkably similar to the MCT analysis. The FTSPD yields 
the following kinetic equation for the evolution of the normalized density autocorrelation function in the case of 
Smoluchowski Dynamics: 

where the memory function kernel is given at second order by 

p4j32i:gB(^,0 = ^y'^[Gpp(^,0G|,°^(|q-k|,0+Gpp(fc,0Gpp(|q-k|,0] • (4.1a) 

In writing down (14. lab we have already assumed that the frequency dependence of the 3-point vertices can, reason- 
ably, be dropped in the long-time limit. This amounts to dropping time derivative terms acting on the density-related 

quantities Gpp, Gpp, Gpp, Gpp where Gpp is the non-interacting correlation function, and the dressed propagators 
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Gpp and Gpp are discussed in some detail in Appendix A. In the following section we demonstrate that this fre- 
quency dependence is irrelevant for our asymptotic analysis. To analyze the asymptotic behavior of this equation 
using the methods of MCT we must examine the Laplace transform of (j){q,t): 

zHq.z) - 1 + ^Hq,z) +Dq^^[pi5^LBB{q,t)]{zHq,z) - 1) = o, 

where we have used the initial condition (j){q,0) = 1 and ^[ZBB{q,t)] is the Laplace transform of I,BB[q,t). Dividing 
by z^{q,z) — 1 we obtain 

I4t^%^ = 1 +Dq'^W^Bsiq,t)] . (4.1b) 
S{q) l-z(j){q,z) 

We see that we obtain essentially the same equation as in MCT but with a different kemel. 
4.1 Long-time behavior of EgBiq^t) 

For the asymptotic analysis of our kinetic equation (I4.1bt . we must determine which terms in the kemel I,gg{q,t) 
are relevant at long times. We begin the analysis with an examination of Gpp, defined in the frequency domain by 

Gpp{q,co) = G'^p^OijGjkOkeGf^ 

-V{q) ^GppGBp Ggp + Gpg Gpp Ggp + Gpg GpgGpp j- , 

where a,/ = V'(^) if / 7^ j and is zero otherwise In the time domain this becomes 

Gpp {q, t) =V{qf f_Js ds' {cj,"; (t s' s)Gbp {s)Gfp (s') + G^l (t - s' - s)Gpp (s)g'§1 (s') 



{t~s'-s)GpB{s)Gfp{s')}, 

where we have suppressed the wave-number dependence on the right-hand side. We shall now examine each of the 
three terms in the integrand in turn. Using the fluctuation dissipation theorem, we have for the first term 

r\[jsfjs'Gfp{t -s'- s)GBp{s)Gfp{s') = j'^Js fjs'Gfp {t-s'- s)ei-s)Gppi-s)e{-s')Gfl{-s') 

= fjs j\s'Gfp{t-s'-s)Gpp{-s)Gfp{-s') 
= ^4p(0) fjsGfp{t-s)Gpp{s), 

where we have performed the integral over s' in the final line. For the second term we have 
£^ ds ds' G% {t-s'- s)Gpp {s)d^l {s') 

= jy^ y'"^//0(f-.'-.)G(''^(f-.'-.)GppW0(-,o4p(-^') 

= fjs'f'^ dsG^^^p{t-s'-s)Gpp{s)G^°l{-s') 

" ds' r'\sGfht-s'-s)Gpp{s)d°l{-s')- r ds'd°h0)Gpp{t-s')G^^h-s') 



ppyi-A - A)\Jpp[^AJ\Jpp\^-A ) - J UA \Jpp(\J)\Jpp\l — A )\Jppy 

-Dq^ ds' y'~' dsGfp {t - s' - s)Gpp {s)Gfp {s') + g|,°^ (0) cj,"^ (f - s)Gpp {s) ds 



Gfp{QfGpp{t). 
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Finally, the third term becomes 

j3-2 £° ds P ds' (f - / - s)GpBis)dpl is') = £ ds £° ds' Oit-s'- s)Gf^ it-s'-s)0 {s)Gpp {s)Gf^p {s') 

= [ ds J'~^ ds' Gf^ {t - s' - s)Gpp {s)Gf^ is') 
= -Dq^ £ ds I' ' ds' Gfl (t~s'- s)Gpp {s)dpi is') . 

Next, using the fact that Gpp (f ) = ±Dq^G^pp (f ) to add up the double integrals from the second and third terms (with 
the sign in front of depending on the sign of t), we obtain 

- V £'d, (jf^ (t-s'- s)Gpp (.)g|,°p' [s') - Dq^ £ds pds' Gfl {t-s'- s)Gpp {s)Gf^ {s') 
= ~Dq^ ds 1^ ' ds'Gf^ {t-s^ s')Gpp {s)Gf^ is') - cj^p' (0) £ cj^p' {t - s)Gpp (s) ds 



\d^X0)l\sG^^^pit~s)Gpp{s). 



2 

Adding all three terms, 

Gpp{q,t) = [pViq)]^l^-Dq^ l^dsj^ ds'd°^p{t - s - s')Gpp{s)Gf^pis') 

-d^^O) J^d^^p{t-s)Gpp{s)ds + d^^p{OfGpp{t 

= [l5Viq)]'d^^p (0) I (O)Gpp {t)-j'^ds[\+ Dq\t - s)]Gfp (f - s)Gpp {s) | . 

From here we compute the Laplace transform, 

Gpp{q,z) = m{q)f-Gfp{t = 0) {cfpit = Q)Gpp{z) ^i±^(zGpp(z) - Gpp{t = 0))| . 

We see that the only singularities as z — > come from Gpp (z). Therefore, in the z — > limit we find that the most 
dominant term is 

Gpp{q,z)~^[^V{q)p^YGpp{q,z), 
where we have replaced Gp*]] (0) with po. Taking the inverse transform, 

Gpp{q,t > to) = [l5V{q)pofGpp{q,t) , 

where fo is a timescale governing the early-time dynamics. Plugging in our ansatz, Gpp{q,t) = S{q){fq + hq\t\"'), 
and multiplying by Gpp (f), 

^[Gp,(«,l)G|,p(«.')l(z) = [PV(,)po]=i(<,)Ji?"[(/, + V°)Gl,p(')) 

Since this term is not singular as z — > it does not affect the asymptotic analysis in the lowest order limit. Thus 
we can ignore the "single particle" contribution to L^b in ( 14. lab . Note that the derivative terms in the kernel add 
multiples of z to our Laplace transformed expressions which implies that they go to zero in the limit. It follows that 
the dominant contribution to the asymptotic equation from the kernel is the "collective" Gpp {q — k,t)Gpp {k, t) term. 
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4.2 Ergodic/Non-ergodic Transition 

To determine if an ENE transition occurs we set (l>{q,z) = fq/z and let z — > 0. Plugging into ( I4.1bb . we obtain 

W)^, = 1 +Dq^np^^^BB{q,t)] 
with ^ ^ 

^W^Lbb] = \-^j ^V{k)S{k)V{q-k)S{q-k)hfq.k- 
Putting this form into the equation, and setting 

j32 f d^k 



2p J (27r) 

we obtain 



Fq[fk]=S{q)j- f T^V{k)Sik)V{q~k)S{q-k)Mg-k 



Zl-fq Dq^^z " 

which becomes 

= F,[fk] 



iq 

in the z — ;> limit. From here it is apparent that the determination of the ENE transition is identical to MCT, except 
that we have a different "vertex" or potential. 

4.3 Time Dependence at Transition 

Ignoring the constant term (which becomes negligible in the z limit), equation ( I4.1bl ) reduces to 



l-z(/>(^,z) 

We apply the ansatz ^{q,t) = f{q) + (1 — f{q)Y'^q{t). The Laplace transform of the ansatz is given by 

(/.(^,z) = ^ + (l-/(^))V,(z). 

z 

We now plug this into the equation and expand the left hand side to second order in ^fq{z): 

YTj^ (7 + (i-/?)v^?W+z(i-/'?)V,'(z)) = np^^^BB{q,t)] . 

To apply the ansatz to the right hand side we must analyze Gpp. Ultimately, we need to compute the Laplace 
transform of Gpp {q, t)Gpp {p, t) . A prerequisite for this computation is Gpp {q, t) . Our ansatz will be 

Gpp{q,t)=S{q) + 

This implies we set ^rq{t) = hq\t\"- /{\ ~fq)^. 



V{q) 



Gpp (q, (0) = ^ [{Gfp {q, ft)) {GpB{q. (O) + Gsp (q, ft))) + Gpp (q, ft)) (cj'j , «) + G^^J (q, ft))) 
In the time domain, we have (suppressing wave-number arguments) 

Gpp (0 = ds [{Gfp is) (Gps it s) + Gbp [t -s))+ Gpp [s) [G^, {t-s)+ G^^l [t - s)) 

= I [ Gpp [s - t)Gfp {s) ds + J'^ Gpp {t - s)Gf^ {s) ds 



P'' l2Gpp{t-s){G^^j,is) + d^lis))ds\ (4.2) 
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where we have invoked the fluctuation dissipation theorem (FDT) in the second line. Using the FDT on one of the 
terms in the third integral we obtain 

r'J^^ Gpp {t - s)Gfl is) ds = Gpp [t - s)Gfp{-s) ds 

= -Gpp {t)G'^pp (0) - Gpp it - s)dp°^pi-s) ds 

= -Gpp it)GPp (0) - l°^Gpp it - s)dp^^p is) ds , 

where we have used the properties of G^ppit) in the second and third lines. A similar computation reveals that the 
other term in the third integral of (14.21 reduces to 

Adding all the pieces together, 

Gpp (q, t) = I Gpp (s - Og|,°p' is) ds + £^ + £- Gpp (t - s)Gf^p (s) ds - 2Gpp (f )Gpp (0) | 

= pViq) I -Gpp{t)Gfl{0) + 1^ Gpp it - s)Gfl is) ds^ . 
Taking the Laplace transform, 

Gpp iq,z) = Pviq) {-Gpp iz)G^p^p (f = 0) + izGpp iz) - Gpp it = 0))Gf^ (z) } . 

Since G^pp iz) is not singular as z — !■ 0, we find that the non-integral term of Gpp it) dominates in the long-time limit. 
Therefore, we may replace 

Gpp(^,0 ^ -l5poViq)Siq)if, + h,\tn 

in this limit. At this point the analysis becomes identical to the MCT case, except that we have a different vertex 
function. Using the notation of Section 3, let 

m] = / ^Vik)Sik)Viq-k)Siq-k)M,., 



5p 

where Viq) is a pseudo-potential determined in 0]. Putting all this together, we find 



limif [i:BB(^,f)](z) = ^^+C,,iifkiz)+Cgkp^[iifkit)Vp{t)] ■ 
z^O Z 



Therefore, 



Fqlfk] + {Wc,iz)-C,kYk{z)) +zil-f,)Y,iz)'C,,p^[Ykit)Ypit)]=0 



and the analysis for the ergodic and non-ergodic regimes now proceeds identically to the mode coupling case, the 
only difference being the "vertex" in MCT and the pseudo-potential in FTSPD. 
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Fig. 7 The product V{k)S{k) appearing in the FTSPD kinetic equation kernel EsBigj) = / ■^^ykSkVci-kSq-k^{k,t)^{q - k,t) at 
the transition. The negative values of this combination imply that C^J is not a positive matrix. This is why we cannot choose such 
that > for all k as we can in Mode Coupling Theory. Nevertheless, the maximum eigenvalue of C^J approaches unity as the packing 
fraction tends to the transition point implying the analysis is still valid. 

5 Numerical Results for Hard Spheres 

Solving for / and iff as suggested by the preceding MCT analysis yields results in close agreement with UbI i. We 
find that with the Percus-Yevick approximation for the structure factor of hard spheres the transition happens at a 
packing fraction of t]mct = 0.515 and Amct = 0.734. Solving for the exponents, we find a = 0.312 and b = 0.584. 
We find that A increases quite sharply in the vicinity of the critical density. Barrat et. al. report t]mct = 0.525 and 
■^MCT ~ 0.758 with the Verlet- Weiss approximation. 

While conducting these numerical experiments it became clear that convergence near the critical density requires 
integrating out to fairly high {ko ^ 100) wavenumbers. Increasing the range beyond this does not have a noticeable 
effect on the calculations. In addition, the number of iterations required for the convergence of the non-ergodicity 
parameter increases quite rapidly as the transition density is approached. For instance, convergence to three digits 
in the eigenvalue occurs after 60 iterations for t]mct = 0.52, but it takes over 700 iterations at 77MCT = 0.51505. 

The vertex in MCT has some special properties which allow one to draw conclusions about the eigenvalues and 
eigenvectors. Because C^''^ is a positive matrix it can be proved that there is a largest eigenvalue and that it becomes 
unity in the continuum limit. In addition, it can be shown that the both the right and left eigenvectors, e^^ and ej., 
must be positive. In FTSPD the situation is more complicated - as can be seen from Figure [Tjthe pseudopotential no 
longer guarantees that C^^ > 0. Therefore the same proofs from MCT cannot be carried over to the FTSPD case. 
However, despite the complication with the pseudopotential, numerical evidence (see Figure |8ll indicates that the 
largest eigenvalue of C^j^ is unity. In addition, the right eigenvector ei^is found to be positive, as in MCT. We also 

observe the same s/e scaling behavior in f^^ that we see in MCT (see Figure |9]l. 

Using the first-order result for the pseudopotential for simplicity, which is essentially the direct correlation func- 
tion fll,'?], solving the FTSPD for hard spheres with the Percus-Yevick approximation appears to be computationally 
easier than MCT. Convergence only requires integration to at most ka ^ 80 compared to ka ^ 100 for MCT. This is 
expected because the decay of the vertex function is much faster for large wavenumbers in the FTSPD. Convergence 
requires more iteration closer to the transition similarly to MCT. We find that the ENE transition occurs at a packing 
fraction of t]ft = 0.602 with a corresponding Aft = 0.553. The non-ergodicity parameter at the transition is plotted 
in Figure [10] while the right and left eigenvectors are plotted in Figures [TT] and [12] respectively. Solving for the 
exponents, we obtain a = 0.386 for the power law relaxation and b = 0.935 for the von Schweidler relaxation. 
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Fig. 8 The eigenvalues of sorted by eigenvalue modulus. The largest eigenvalue is real, non-degenerate, and approaches unity as 
the density approaches the transition point. This graph was constructed from C^J computed on a 800 x 800 grid with a spacing dk = 0.1. 
The spectra for grids down to 400 x 400 with dk = 0.2 and up to 1600 x 1600 with dk = 0.05 are essentially identical. 



6 Conclusion 

We have demonstrated that, to second order in the effective potential, the FTSPD essentially reduces to a mode 
coupling theory. This shows that FTSPD is a framework from which MCT can be derived and systematic corrections 
introduced. We find that an ergodic / non-ergodic transition occurs at a packing fraction of 77ft = 0.602, which is 
in good agreement with some recent simulations and experiments [14]. This agreement prompts us to investigate 
the effectiveness of the pseudopotential approach for other systems. In addition, we plan to confirm this asymptotic 
analysis by a direct integration of the FTSPD kinetic equation. Although we have gathered extensive numerical 
evidence for their veracity, proofs that the largest eigenvalue Eq of is real, positive, and nondegenerate and 
£o — ^ 1 as 77 — ;> 77 * are still unavailable to the best of our knowledge. In a forthcoming paper we extend this analysis 
to a treatment of binary mixtures. An analysis of Newtonian dynamics within the framework of FTSPD is given in 

m. 



A Properties of G and G 

In this appendix we prove that the dressed propagators G and G individually satisfy the same fluctuation dissipation 
theorem that G satisfies: ^ 

Gpp {q, ®) = ^ Im Gsp {q, ft)) . 

In the time domain, 

GBp{q,t)=Pd{-t)Gpp{q,-t). 
This result holds at all orders of perturbation theory. The non-interacting cumulants are given by 



(0). IpoDq 



2 



G';>iq,CO) = 



and 



Gp,iq,CO)-G,p{q,CO) - ■ 
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Fig. 9 This plot shows f^{ka = 7.4), the non-ergodicity parameter at wavenumber ka = 7.4, as a function of the packing fraction. 
The dashed line shows a least squares fit of A^/e+B to the data. 
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Fig. 10 This plot shows fF at the critical density ?Jft = 0.602. 



In the following subsections it will be useful to decompose the cumulants into real and imaginary components: 



d^l^Ro + ilo 
gS=/?o-«7o 



■ R + il 



GpB =R — iI 



15(0 
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Fig. 11 This plot shows ej^, the right eigenvector of C^, at the critical density rjpT = 0.602. The eigenvector is normaUzed according 
to the convention J^i^ekiii = 1. 

A.l GFDS 

Recall the form of given by 



Explicitly, this yields 



1 



^(0)x_^(0), 



Gbp = l{Gfll5VGBp + GBpl3VGfl) = G^bI^VGbp , 



and 



Gpp = UciyVGpp+GflliVGBp + GpBliVGf^+GppliVGfl). 



If we write out our contributing terms as real and imaginary components then we have for the imaginary part of 

GpB, 

IrnGsp = Yi^^^^ ~GpB) 



Looking next at Gpp, we have 



= ^ [(^0 + ik) {R + il) - {Ro - ik) {R - il)] 
= pV{IRo+IoR). 



2V 



Gpp = —iRoI+RIo). 



Finally, we obtain the fluctuation dissipation relation 



Gpp = J^^^Bp- 
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Fig. 12 This plot shows ef^, the left eigenvector of C^J, at the critical density rjpx = 0.602. We see that this eigenvector is not positive for 
all k, in contrast to MCT. The reason for this is that is not a positive matrix, rendering the Perron-Frobenius theorems inappMcable. 
The eigenvector is normalized according to the convention Y,k^kek = ^- 



A.2 G FDS 

We may repeat the same procedure for Gij, which is given by 
This implies 



GpB = Glp^d^ll5VGpBPVG^^j, 
and 

Gpp = GPpPVGBpPVG^^l+G^^^sPVGppPVG^^l+GflpVGpBPVG^^^ . 
Using the same decomposition into real and imaginary components, we have 

Im Gbp = — {Gbp — GpB) 

r; (G^Bl^yGnp^yG^^l-G%^VGpB^VG% 



2i 



[{Ro + Ho) {R + il) {Ro + Ho) - {Ro - Ho) {R - H) {Ro - Ho)] 



and 



2i 

mf{{Rl-ll)I + 2RoIoR) 



Gpp = [lo {R + H) {Ro + Ho) + {Ro + Ho)I{Ro - Ho) + {Ro - Ho) {R - H)Io] 



[{Rl-ll)I + 2RoIoR] . 



(O 



Therefore, 



2 
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B Conventions and Fundamental Relationships 

We define the Fourier tranform of a function f{t) as 

J — oo 

with inverse transform 

One of the tools we have used extensively is the Laplace transform: 

^[f{tm= re-^f{t)dt. 

Jo 

We define the convolution of two functions / and g as: 

{f*g)it)= f f{t-s)g{s)ds 
Jo 

The Laplace transform of a convolution is then 

^[if*8m]=j^[f{t)]^[8{t)]. 

In addition we have employed a "discrete" notation for various integral kernels. This is convenient since we are 
only considering spherically symmetric functions for which the notation corresponds to the dicretized form used in 
the numerical routines. For example, 




and 

5,fc = (2;r)3 5(q-k). 
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